Topological phases in gated bilayer graphene: Effects of Rashba spin-orbit coupling 

and exchange field 
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We present a systematic study on the influence of Rashba spin-orbit coupling, interlayer potential 
difference and exchange field on the topological properties of bilayer graphene. In the presence of 
only Rashba spin-orbit coupling and interlayer potential difference, the band gap opening due to 
broken out-of-plane inversion symmetry offers new possibilities of realizing tunable topological phase 
transitions by varying an external gate voltage. We find a two-dimensional Z2 topological insulator 
phase and a quantum valley Hall phase in AB-stacked bilayer graphene and obtain their effective 
low-energy Hamiltonians near the Dirac points. For AA stacking, we do not find any topological 
insulator phase in the presence of large Rashba spin-orbit coupling. When the exchange field is 
also turned on, the bilayer system exhibits a rich variety of topological phases including a quantum 
anomalous Hall phase, and we obtain the phase diagram as a function of the Rashba spin-orbit 
coupling, interlayer potential difference, and exchange field. 

PACS numbers: 73.22.Pr, 73.43.Cd, 71.70.Ej, 73.43.-f 



I. INTRODUCTION 

Topological insulator— (TI) is a new phase of quan- 
tum matter in materials with strong spin-orbit cou- 
pling. The quantum spin Hall effects in graphene 2 - and 
HgTe quantum wells 3 - represent the first examples of two- 
dimensional topological insulators. So far, semiconductor 
heterostructures like HgTe^ and InAs/GaSb 6 quantum 
wells offer the only realistic materials with strong spin- 
orbit coupling that can realize the quantum spin Hall 
phase, and the intrinsic spin-orbit coupling in graphene 
was shown to be too weak. — Due to graphene's attrac- 
tiveness as a potential electronic material for emerging 
nanotechnology, artificially enhancing the spin-orbit cou- 
pling strength in graphene can open up new possibili- 
ties in graphene-based spintronics, and several theoreti- 
cal and experimental works have addressed the effects of 
enhanced spin-orbit coupling in graphene by doping with 
heavy adatoms such as indium or thallium 7 , doping with 
3d/5d transition metal atoms^— , and interfacing with 
metal substrates, e.g., Ni(lll)ii. 

Due to band gap opening from broken out-of-plane in- 
version symmetry, gated bilayer graphene is a quantum 
valley-Hall insulator (QVHI) characterized by a quan- 
tized valley Chern number. In our recent work^ we 
have reported that the presence of Rashba spin-orbit cou- 
pling turns the gated bilayer graphene system from a 
QVHI into a Z2 TI, with the phase boundary given by 
X 2 R = U 2 + 1\ where Xr, U and t± denote the strengths 
of Rashba spin-orbit coupling, interlayer potential differ- 
ence and interlayer tunneling amplitude respectively. In 
this paper, we obtain low-energy effective Hamiltonians 
for the topological insulator phase, valid for small U <C 
Xr and below the topological phase transition, as well as 
for the quantum valley Hall phase above the phase tran- 
sition. In the presence of different Rashba spin-orbit cou- 



pling strengths on the top and bottom layers (X R ^ X R ), 
we show that the topological insulator phase remains ro- 
bust as long as X R X R > t\. When the time-reversal sym- 
metry is broken by an additional exchange field M, the 
bilayer system hosts different topological phases charac- 
terized by different Chern numbers C = 2,4sgn(M) and 
valley Chern numbers C v — 2,4sgn(Z7), and the phase 
boundaries associated with the topological phase transi- 
tions are given by U = ±M and U 2 + 1\ - M 2 - X 2 R = 0. 

The rest of this paper is organized as follows. Sec- 
tion [II] introduces the tight-binding and low-energy effec- 
tive Hamiltonian of Ai?-stacked bilayer graphene in the 
presence of Rashba spin-orbit coupling, exchange field, 
and interlayer potential difference. In Section IIIII we ob- 
tain the low-energy effective Hamiltonians of the Z2 TI 
and the quantum valley Hall phases. We consider in Sec- 
tion IIVI the case of AA-stacked bilayer graphene. The 
effect of different Rashba spin-orbit coupling strengths 
on the top and bottom layers of ^4i?-stacked bilayer 
graphene is considered in Section [V] Finally, in Sec- 
tion IVII we obtain the phase diagram of AB-stacked bi- 
layer graphene as a function of X R , U and M. 



II. SYSTEM HAMILTONIAN 

The tight-binding Hamiltonian of A£>-stacked bilayer 
graphene in the presence of Rashba spin-orbit coupling, 
exchange field, and interlayer potential difference is given 
byi 2 .^ 3 



#BLG - #SXG + H SLG 
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where the first two terms H^q denote the monolayer 
graphene Hamiltonian for the top (T) and bottom (B) 
layer (see below) , the third term represents the tunneling 
Hamiltonian that couples the top and bottom layers, and 
the last two terms take into account a potential difference 
2U between the top and bottom layers. The single-layer 
graphene Hamiltonian is& & 14 i 15 



with 



Hslg — H o + Hr + Hm , (2) 

(H) 

H R = it R ^ e z -(s a 0Xdij)cl a Cjp; 
H M = M c \ a sZ a p Ci ^ 



where (...) runs over all the nearest neighbor sites with 
hopping amplitude t = 2.6 eV, s are spin Pauli matrices 
with a and j3 denoting up spin or down spin, and cj a 
(ci a ) is the electron creation (annihilation) operator on 
site i. Hr describes the Rashba spin-orbit coupling with 
coupling strength tn and is a lattice vector pointing 
from sites j 'to i, and Hm the exchange field contribution 
with magnetization M. 

In the momentum space—, we perform an expansion of 
the momentum about the valley points K and K' and ob- 
tain the following eight-band effective Hamiltonia m 12 ' 13 

H =v{qO X k X + <7yky)l S l T + -^-(<T X T X - <JyTy)l S 

+ ~^-(WxSy ~ <TyS x )l T + MS Z l a l T + UT z l s l a , (3) 

where r\ = ±1 label the valley degrees of freedom, <r and 
r are Pauli matrices representing the AB sublattice and 
top-bottom layer degrees of freedom, 1 is a 2 x 2 identity 
matrix. The bare graphene Fermi velocity is given by 
v = Sat/2 with a the lattice constant and Rashba spin- 
orbit coupling is given by Xr = Mr. For simplicity, we 
set the lattice constant a to be unity henceforth. 



III. FOUR-BAND EFFECTIVE HAMILTONIAN 

In Ref. [12] ]. we reported numerical tight-binding cal- 
culations showing that ^4_B-stacked bilayer graphene un- 
der external interlayer potential undergoes a topological 
phase transition from a QVHI to a two-dimensional Z2 
TI. Figure [1] shows that, at the phase transition point, the 
bulk band gap of the bilayer graphene system obtained 
from Eq. ([3]) closes at the valley points at the zero energy 
(for comparison Figure []] also shows the band structure 
of a pristine single-layer graphene). It is therefore pos- 
sible to obtain low-energy Hamiltonian descriptions near 
band gap closing, which occurs when bilayer graphene 
turns into a TI from a semi-metal for small U and at the 




FIG. 1. (Color online) Solid line: bulk band structure of 
bilayer graphene at the phase transition point between the 
quantum valley-Hall and the Z2 topological insulator phases. 
The parameters used are U = O.lOt, t± = 0.143t and t_R = 
0.058£. Dashed line: bulk band structure of pristine single 
layer graphene. The band gap closing at phase transition 
occurs precisely at the valleys K and K' . 
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FIG. 2. (Color online) Phase diagram of bilayer graphene 
as a function of U, Xr, t±. The region inside the ellipsoid 
corresponds to the Z2 topological insulator state and outside 
region to the quantum valley Hall insulator phase. The phase 
boundary is given by Eq. Q 



phase transition point U = Uq between TI and QVHI. 
In the following, we expand the eight-band Hamiltonian 
Eq. ([3]) up to the leading order in momentum k around 
the valley points and obtain a reduced four-band effec- 
tive Hamiltonian that captures the low-energy physics of 
the system near phase transitions. 

We assume equal Rashba spin-orbit coupling strengths 
in both top and bottom layers in this Section and address 
the effects of unequal Rashba strengths in Section[V] Fig- 
ured] illustrates our results for the energy bands obtained 
numerically from the eight-band Hamiltonian in Eq. ([3]) . 
We observe that the fourth band and the fifth band be- 
come inverted when the Rashba strength increases be- 
yond a critical value signaling a topological phase tran- 
sition. Imposing k = in Eq. (PJ gives the condition for 
gap closing of the bulk bands as 12 

\l = U 2 + tl. (4) 

The (U, Xr, t±) phase space is therefore divided into 
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FIG. 3. Energy dispersion at valley A" (k = 0) as a function 
of Rashba spin-orbit coupling Ar at fixed U = 0.3t and t± — 
0.143i. One notices that band inversion occurs between band 
'4' and band '5' at the critical Rashba spin-orbit coupling 
value Xr ~ 0.33£. 



two regions as illustrated in Fig. [5] the system is in the 
C v = 4sgn(C7) QVHI phase when X 2 R <U 2 + t\, while it 
is in the Z 2 TI phase characterized by C v — 2 sgn({7) and 
Z 2 = 1 when A| > ?7 2 + 

In Fig. 01 we plot the Chern number contributions from 
each valence band using the eight-band Hamiltonian in 
Eq. © near the valley K. Here, we fix the Rashba spin- 
orbit coupling and interlayer coupling as Xnjt = 0.617 
and t±/t — 0.143. The topological phase transition point 
occurs at Uq 



t,. In the TI phase when U < Uq, 



the plot shows that the contribution to the total Chern 
number from each valence band varies as a function of U, 
and in particular there are two regimes where the Chern 
number contributions are distributed differently among 
the bands. For U 0, Cf 2 3 4 = -0.5,0.5,-1.5,0.5 

while for U -> C/ " , Cff 2 ,3,4 = 0.0, 1.0, -2.0,0.0. For the 
QVHI phase occurring when U > Uq, the Chern num- 
ber contribution from each valence band is constant as a 
function of U with Cf 2 3A = 0.0, 1.0, -2.0, -1.0. In the 
following we study the low-energy physics of the TI and 
QVHI states in these three regimes. 



Near Semimetal-TI Phase Boundary: U — > and 

U<X R 



FIG. 4. (Color online) Chern number contribution from each 
valence band obtained from Eq. [3] at valley K as a function 
of U /Uo. The Rashba spin-orbit coupling and interlayer tun- 
neling are fixed. Co is the critical value (see the dashed line) 
of interlayer potential separating the 2D TI and QVH phases. 
For U < Uo, the total Chern number Ck = — 1, but the 
contribution from each valence band varies as a function of 
U /Uo- For U > Uq, the total Chern number Ck = —2 and the 
contribution from each valence band is independent of U /Uq. 
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In the limit U — Y and U -C Xr, Hi and H2 corre- 
spond respectively to the lower bands [i.e. e = ±E7] and 
higher bands [i.e. e = ±(Ar ± U)]. In the vicinity of K, 
the coupling T between H\ and H% becomes very weak. 
Therefore, the original eight-band Hamiltonian can be 
reduced to an effective four-band Hamiltonian^ 



In the basis {A lt , B n , A 2 f, B 2 i, A n , B lt , A 2 i, B 2 f}, 
the eight-band Hamiltonian Eq. ([3]) can be written at 
valley K as: 



H K = 



H l T 
T H 2 



(5) 
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FIG. 5. Bulk band structure of the effective four-band model 
in the limit of Xr/U 3> 1 along the profile of k y = 0. '1' and 
'2' label the two valence bands. Here, we choose Xr/ = 0.20t 
and U — 0.02t. Note that the two conduction or valence 
bands touch at k = 0. 



Similarly, the eight-band Hamiltonian at valley K' in 
the basis of {Ai^ 
can be expressed as 



-Bit' -^24.) B 2 \, Ai^, B11, A 2 f, B 2 i} 
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Using the formula of Eq. ©, the resulting reduced 
four-band Hamiltonian can be written as: 
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Upon diagonalization of Hamiltonians in Eqs. © and 
(fT2"j). the energy dispersions at both K and K' can be 
obtained and share the same form: 



e = ± 



yJ\* R U* + e\ + 2e 2 t\ ± 2^P^TI 



2j-2 
fc fc _L 



(13) 



where — vk. 

Figure [5] plots the bulk band structure of the four-band 
effective Hamiltonian in Eq. (|13|) along k v = 0. One can 
see that a bulk band gap opens and the resulting conduc- 
tion and valence bands touch at k x =0. The bulk band 
gap opening signals an insulating state. To reveal its 
topological property, we have evaluated the Berry phase 
contribution from the occupied valence bands below the 
band gap. In a continuum model Hamiltonian, the Chern 



number is calculated by integrating the Berry curvature 
in the entire momentum space: 



1 
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dk x dky^l n {k x , ky)^ (14) 



where VL n {k x ,k y ) is the momentum-space Berry curva- 
ture at (k x , k y ) of the n-th band, and is given by 



E 2Im( 1p nk \v x \l/> n > k ) (4>n> k\Vy\lpnk) , .. 
77. 77 — ' \ Lb > 



where u n = £ r JTi, and v x t y \ is the velocity operator along 
the x(y)-direction. 

In Fig. [6j we display the Berry curvature distribution 
f2 along k y = for both valleys K and K'. One ob- 
serves that the Berry curvatures are exactly opposite at 
the two valleys K and K' for each band. In particular, we 
find that the total Berry curvatures around K or K' do 
not share the same sign in the whole momentum space, 
in contrast with the Berry curvature distribution in the 
quantum anomalous Hall effect in single layer graphene 9 
or the conventional QVHI in graphene due to the pres- 
ence of staggered AB sublattice potential . 10 ' 19 

By numerically evaluating the integration in the mo- 
mentum space, the Chern numbers at valleys K and K' 
are found to be 



C K = -C K > = sgn([/), 



(16) 



in which the Chern number contribution from each va- 
lence band arc 



r 1 



r 2 



-C 2 



K' 



;Sgn(i7), 



(17) 
(18) 



where the supscripts label the valence band indices in 
Fig. El 

From the principle of bulk-edge correspondence, one 
expect that there is only one pair of edge states propa- 
gating along the system boundaries. Since time-reversal 
symmetry is preserved in the our system, one concludes 
that this nontrivial insulating state belongs to the Zi 
TI class. Therefore, the edge modes are robust against 
weak non-magnetic impurities. Moreover, since differ- 
ent valleys are encoded into the counter-propagating edge 
channels^ they are further protected by the large mo- 
mentum separation as long as inter-valley scattering is 
forbidden. As a consequence, these edge modes are also 
robust against smooth non-magnetic and magnetic im- 
purities. From Eq. (|16p . the valley Chern number of this 
Z 2 TI state is: 



2e 2 

C V =C K - C K > = — sgn(f7), 
h 

which is half of that in the conventional QVHI 
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FIG. 6. Berry curvature distribution Q at valleys K and K' 
along k y = of the effective four-band model in the limit of 
Xr/U 3> 1. Upper panels: Berry curvature distribution Qi 
for the lowest valence band labeled as 'I' in Fig. Middle 
panels: Berry curvature distribution fia for the valence band 
close to the band gap labeled as '2' in Fig. [S] Lower panels: 
the total valence band Berry curvature distribution fir- The 
parameters adopted here are the same as those used in Fig. [5] 



B. Near TI-QVHI Phase Boundary: U ->■ ^JX 2 R - t]_ 

In the following discussion, we set Xr and t± as fixed, 
and allow U to slightly deviate from the p hase tran sition 
point U , i.e., U = U + A, where U = \/A|j — A > 
and < correspond respectively to the QVHI and TI 
phases. 

When the topological phase transition occurs, the bulk 
band gap closes at the valley points K/K' . Hence, at the 
critical value Uq and k = 0, the system Hamiltonian in 
the basis of {B n , A 2 f, A it , A 2i , B 2 ^, B 2i , B lt , A n } can 
be expressed as 
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The eigenenergies can be obtained from the above as e — 
0, 0, ±U , ±e 1: ±e 2 , where £l = -Z7 /2 + ^/8A^ + 
and e 2 = -U /2- y/8X 2 R + Ufi/2. The former four cor- 
respond to the low-energy part, while the latter four cor- 
respond to the high-energy part. Based on our analy- 
sis of the Berry curvature from the tight-binding model, 
the valley Chern number arises only from the low-energy 
bands and the high-energy bands give no contribution. 

The unitary transformation matrix that diagonalizes 
Eq. (|2"0")) is presented in the Appendix. After some ma- 
nipulations, the effective Hamiltonian can be obtained 



H cS = H P -TH Q 1 T\ (21) 

where explicit expressions of Hp , T, Hq are also given in 
the Appendix. The effective Hamiltonian to first order 
in k is then given by 



H. 
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V2X R 

In the following, we show that the first-order effec- 
tive Hamiltonian at valley K is sufficient to capture the 
topological phase transition between the QVHI and TI 
phases. Figure [7] displays the band structures along k y = 
for three different values of A/t = -0.03,0.00,0.03 at 
fixed Uo/t = 0.30 and Xp/t = 0.33. One observes that 
when A increases from negative to positive, closing and 
reopening of the bulk band gap occur as expected. Fig- 
ure [5] plots the Berry curvatures along k y = of the two 
valence bands below the band gap and the total Berry 
curvature for A/t = ±0.03. The Berry curvatures for 
the first valence band fii in both cases are similar shar- 
ing the same sign [panels (a) and (d)]. For the second 
valence band, the Berry curvatures for A/t = —0.03 ex- 
hibit both positive and negative signs [panel (b)], but for 
A = +0.03, they are both negative [panel (e)]. As a con- 
sequence, the total berry curvatures at A/t = —0.03 in- 
clude both positive and negative contributions [panel (c)], 
and those at A = +0.03 share the same negative sign. By 
using Eq. (|14[) . the Chern numbers acquired by the va- 
lence bands at valley K is evaluated as Ck = —1,-2 
for A = —0.03, +0.03 respectively. Following a similar 
derivation, one obtains the effective Hamiltonian at val- 
ley K' with the corresponding Chern numbers Ck' =1,2 
for A = —0.03, +0.03 respectively. In this way, the val- 
ley Chern number in the topological insulator phase is 
C v = Ck — Ck> = —2 while that in the quantum valley 
Hall phase is C v — Ck — Ck 1 = —4. 

In all the three limits, the valley-Chern numbers from 
the resulting four-band effective Hamiltonians are consis- 
tent with those from the direct eight-band full Hamilto- 
nian. 



IV. AA-STACKED BILAYER GRAPHENE 

Bilayer graphene is composed of two monolayers of 
graphene, usually arranged in AB or AA stacking pat- 
tern. In previous Sections, we have predicted a Z 2 TI 
phase in the AB stacking configuration. It therefore be- 
comes a natural question to ask whether the AA stacking 
configuration can also host a Z 2 TI phase. In this Section, 
we demonstrate that the AA-st&cked bilayer graphene 
does not realize a ^topological insulator state in the 
presence of Rashba spin-orbit coupling and interlayer po- 
tential. 
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A/t=-0.03 A/t=0.00 M=+0.03 




FIG. 7. Band structure obtained from the effective Hamilto- 
nian Eq. (1221) along k y = 0. The parameters for the phase 
transition point are set to be Uo/t = 0.30 and X/t = 0.33. (a) 
A/t = -0.03; (b) A/t = 0.00; (c) A/t = +0.03. The bulk 
band gap decreases toward zero when A/t = 0.00 and then 
reopens when A becomes positive. 




FIG. 9. (Color online) Bulk band structures of AA-stacked 
bilayer graphene along k y = at valley K. (a) Pristine 
graphene, U/t = and X R /t = 0; (b) U = O.Wt and X R = 0; 
(c) U = and X R = 0.15t; (d) U = O.lOt and X R = 0.15*. 
Solid and dashed bands correspond to the top and bottom 
layers, respectively. Here the interlayer coupling is set to be 
t± = O.Wt. 



QVHI 




(b) 



(c) 




-(f) 
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Figure [5] plots the bulk band structures of AA-stacked 
bilayer graphene along k y — obtained from the low- 
energy continuum model Hamiltonian at valley K ^ For 
pristine AA-stacked graphene, the linear Dirac dispersion 
near valley K still holds as shown in panel (a), resem- 
bling two copies of monolayer graphene with a relative 
shift of 2t± (solid and dashed lines are used to label the 
bands from top and bottom layers). In the presence of 
an interlayer potential difference, there is no bulk band 
gap opening [see panel (b)] since the inversion symme- 
try with respect to the graphene plane is not violated. 
If only the Rashba spin-orbit coupling is turned on, one 
finds that again the resulting band structure is a com- 
bination of two copies of the monolayer graphene's band 
structures with a relative shift [see panel (c)]. When both 
Rashba spin-orbit coupling and interlayer potential dif- 
ference are present, no bulk band gap appears. We there- 
fore conclude that inversion symmetry breaking is a nec- 
essary requirement for the Z 2 TI in the bilayer graphene 
system. In addition to AA and BB stacking, twisted 
bilayer graphene presents another possibility which has 
attracted much recent interest In future works, it 
will be interesting to study the possibility of inducing a 
Zi TI in a twisted bilayer graphene. 
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FIG. 8. Berry curvature distribution f2 along k y = for A = 
±0.03. Other parameters are Uo/t = 0.30 and X/t = 0.33. 
(a)-(c): Berry curvature distribution for each valence bands 
and the summary of all valence bands at A/t = —0.03. (d)-(f ) 
Berry curvature distribution for each valence bands and the 
summary of all valence bands at A/t = +0.03. 



V. EFFECTS OF DIFFERENT RASHBA 
SPIN-ORBIT COUPLINGS ON TWO LAYERS 

We have studied the Z 2 TI state while assuming the 
same Rashba spin-orbit coupling in the top and bottom 
layers of the AS-stacked bilayer graphene. In bilayer 
graphene, Rashba spin-orbit coupling is extremely weak, 
and one has to employ external means to enhance the 
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Rashba spin-orbit coupling, e.g., doping it with heavy 
metal atoms or interfacing it with substrates. Therefore, 
the resulting Rashba spin-orbit couplings are likely to be 
different on the top and bottom layers. In the following, 
we discuss the effect of different top and bottom Rashba 
spin-orbit coupling strengths on the resulting TI state. 

We adopt the low-energy continuum Hamiltonian 
Eq. and assume different values of top and bottom 
Rashba spin-orbit couplings Ai ^ A2. We first consider 

I 



the case when both |Ai — A2 and U are small. In this 
case, we found from our numerical calculations that al- 
though the conduction and valance bands are no longer 
symmetric with respect to e = 0, the bands still close 
exactly at the valley points K and K'; it is thus possi- 
ble to obtain an analytic formula that describes the band 
closing condition. After imposing k = 0, one finds that 
two of the eigenenergies are e = ±J7, while the remaining 
six ones satisfy the following equations: 



e 3 + {-iyUe 2 - (A? + U 2 + t\_)e - {-\) l U{U 2 - A? + 1\) = 0,i = 1, 2. 



(22) 




FIG. 10. (Color online) The interlayer potential U that closes 
the band gap at the valley K point as a function of Rashba 
spin-orbit couplings Ai and A2, as predicted by Eq. (I26[) . Col- 
ors represent the amplitude of U. It can be seen that while 
such contours are hyperbolas for small U (see the dotted line 
in the lower part of the graph), they become linear for larger 
U (see the straight lines in the upper-right corner). In addi- 
tion, the constraint of A1A2 > t± is clearly demonstrated in 
this plot. 



We search for the condition when the top of the valence 
band and the bottom of the conduction band touch, clos- 
ing the bulk band gap at the K or K' point. This can 
be translated into the condition that the two equations 
in Eq. (|22|) have a common real- valued solution that lies 
between — U and U. The latter condition is necessary 
in order to rule out the scenario that two higher (lower) 
bands touch at the K point. 

It turns out that the numerical search for a common so- 
lution is not as simple as the case with identical Rashba 
effects in Ref. [IH, where we can directly require the 
common solution to be e = 0. In the present case, how- 
ever, the band closing point is no longer fixed at e = 0, 
which makes it difficult to obtain a simple analytical so- 
lution. Instead, we opt to solve for A.; (and not for e) 



I 

from Eq. (f2"2")) . obtaining 

U ~ £Q (TT 2 

U + e 
, _ U + eo fT - , 
U-e 



a: 



(U'+ti-e 2 ), 
(U 2 +t 2 ± -e 2 ), 



(23) 
(24) 



where £0 is the common solution of the two equations in 
Eq. (|2"2")l . Then the band-closing point can be analytically 
obtained by dividing Eq. ([23]) by Eq. ([24]) : 



£0 



A2 — Ai 
A2 + Ai 



U. 



(25) 



In order to have a band gap closing, these parameters 
must also satisfy the following condition 



i + 



4AiA 2 



(Ai + A 2 y 



-U = AiA 2 



(26) 



which is derived by multiplying Eqs. (|23|) and (|24|) . It is 
reassuring to see that when Ai = A 2 , this condition does 
reduce to the one given in Eq. (|4|). One can also rewrite 
Eq. (|26p by expressing the interlayer potential difference 
U as a function of Ai and A 2 : 



U 



t 2 



A1A2 



(27) 



where we see that the band gap at the K or K' point will 
not be able to close if A1A2 < t\. 

Figure [TU] plots the interlayer potential difference U in 
the (Ai, A 2 ) plane that satisfies the band gap closing con- 
dition Eq. (l2"rJl) . Colors represent the strength of U . In 
the blank region, no band gap closing occurs under the 
constraint of AiA 2 > t 2 ^. In the limit of small poten- 
tial difference U, the contours of U behave as hyperbolas 
given by AiA 2 = t\ (see the black dotted line), while 
in the large U limit, the contours tend to straight lines 
given by Ai + A 2 = 2U . 

To verify the correctness of the analytical expression 
of the band gap closing in Eq. (|26|) , we compare it with 
band gap results from direct numerical diagnolization of 
the eight-band continuum model Hamiltonian at valley 
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FIG. 11. (Color online) Numerically obtained band gap at 
valley K point as a function of the interlayer potential U and 
bottom-layer Rashba spin-orbit coupling A2 at a fixed top- 
layer Rashba spin-orbit coupling \\/t — 0.2. Color measures 
of the size of the band gap. The white dots are plotted from 
Eq. (|26[) . The analytical condition Eq. (|26|) for band gap 
closing overlaps with the numerical result. 




FIG. 12. (Color online) Numerically computed band gap at 
valley K point as a function of Rashba spin-orbit couplings 
in both top and bottom layers Ai and A2. In this plot the 
interlayer potential difference is set as U = 0.20t. Color mea- 
sures the size of the band gap in units of t. The white dots 
are plotted from Eq. (|26|) . The analytical result for the band 
gap closing condition agrees well with that obtained from nu- 
merical calculations. 

K/K' . In Fig. [TTJ we plot the band gap at K point as a 
function of interlayer potential difference U and Rashba 
spin-orbit coupling A2 at a fixed Ai = 0.20£. White dots 
plot the analytic result Eq. (|26l) corresponding to the 
boundary for band gap closing. We find that it agrees 
well with the numerically obtained condition for band 
gap closing. Similarly, in Fig. [T^l we plot the band gap 



at K point as a function of the two different Rashba 
spin-orbit couplings Ai and A2 at a fixed interlayer po- 
tential difference U — 0.2t. The white dots obtained 
from Eq. (f2"rJ|) fit exactly where the band gap closes at 
the valley K/K' point. Therefore, the analytical phase 
boundary in Eq. (126p indeed captures the band gap clos- 
ing condition at valley K/K'. 

For very different Ai and A2 and a large U, we find 
that the conduction and valence bands can close indi- 
rectly at different momenta, and as a result there is no 
global bulk gap even though the direct gap at the valley 
points is nonzero. The global bulk gap is the smallest 
energy difference between the conduction band and the 
valence band across the entire Brillouin zone. 

In Fig. 1131 we show the comparison between the nu- 
merically computed band gap (circle or triangle) and di- 
rect band gap at the valley points given by Eq. (f26|) (solid 
line). One observes that for large differences in Ai and A2 
and for large U, the numerically computed gap deviates 
from Eq. (f2l)|) , indicating that band gaps close indirectly 
at different momenta. 

To examine the nontrivial topology of the phases be- 
fore and after the band gap closing, we calculate the 
valley-Chern numbers using the continuum model and 
the Z2 topological number using the tight-binding model 
presented in Ref. (2(| . As depicted in Fig. [T31 before the 
phase transition, the system hosts a QVHI phase with 
C v = Ck — Ck> = 4, while after the phase transition, it 
enters into a TI phase with Z% = 1. Note that the TI 
state is simultaneously a QVHI state characterized by 
C v = 2. These results are consistent with our findings in 
the presence of identical Rashba spin-orbit couplings. 12 
Therefore, we have shown that the Z2 TI state we pre- 
dicted in Ref. [l2| remains robust when the top and bot- 
tom layer Rashba spin-orbit coupling strengths become 
different. 



VI. EXCHANGE FIELD EFFECT 

This section is dedicated to investigate the exchange 
field effect on the QVHI state and Z2 TI state in gated 
bilayer graphene. In Ref. @, we have found that a bulk 
band gap will open in monolayer graphene in the pres- 
ence of both Rashba spin-orbit coupling and exchange 
field, inducing a quantum anomalous-Hall phase^Zr— . In 
Ref. (HI, we have shown that in gated bilayer graphene, 
when the exchange field M is larger than the interlayer 
potential difference U, i.e. M > U, the system undergoes 
a topological phase transition from the C v — 4 QVHI 
phase into a C = 4 quantum anomalous-Hall phase. In 
the following, we supplement this result with a new phase 
boundary that has been overlooked in Ref. [13( . We also 
discuss how the QVHI state with C v = 1,2 evolves in the 
presence of the exchange field. 

For simplicity, we set the Rashba spin-orbit couplings 
to be the same in both layers in this discussion. We start 
from the low energy continuum model Eq. (©), and con- 
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FIG. 13. (Color online) Comparison between the conditions 
for global bulk gap closing obtained numerically (circles and 
triangles) and local bulk gap closing predicted from Eq. (|26[) 
(solid line), for interlayer potential difference U = 0.05t and 

u = o.iot 



sider the following ingredients in our system: interlayer 
potential difference U, Rashba SOC in both layers Xu, 
and exchange field M. The conduction and valence bands 
are symmetric about e — 0, and the bulk band gap closes 
at exactly the valley K/K' point. The energy dispersion 
of the eight-band Hamiltonian at valley K is determined 
by the following equations: 



e = ±(M - 
X 2 R (M±e 



to. 

-U) 



{M±e + U)[U 2 + t\-{M- ef 



The equations for valley K' can be obtained from the 
above by replacing C/ — >- — C/. By imposing e — 0, we 
obtain the following bulk gap closing condition 



M 



U 2 = 0, 



(28) 



which has been reported in Ref . |l3| . In addition, if Xr ^ 



±t± and U 2 
condition 



ti 



A 2 * 



U 2 + t 2 , 



M 



\ R = 0, 



(29) 



is also possible. These two conditions signify two topolog- 
ical phase transitions and give the corresponding phase 
transition boundaries. 

In the above discussions, we have omitted a very im- 
portant case in the presence of vanishing Rashba SOC 
Xr = 0. It is known that a large exchange field will close 
the bulk band gap induced from potential difference and 
results in a metallic phase. Most importantly, the band 
gap closing is not always exactly at the valley points as 
in other phase transitions we have discussed. After rear- 
ranging the Hamiltonian, the 8x8 Hamiltonian can be 
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FIG. 14. (Color online) Phase diagram in the (M, Ah, U) 
space. The whole phase space is divided into various topo- 
logical phases. Here, we only show the phase boundaries. (1) 
The mutual vertical planes determined by Eq. (1281) separate 
the quantum valley Hall and quantum anomalous Hall phases. 
The uniparted hyperboloid determined by Eq. (|29[) is used to 
separate different sub-phases. Note: in the plane of A_r = 0, 
the red color represents the metallic phase. 



written as: 



where 



H 



K 



H x T 
T H 2 



U + M 

U-M 





H 2 = 



U-M 








U + M 








-U + M 





-U-M 








-U - M 







-U + M 



(30) 



(31) 



,(32) 



> then a second gap closing an( j 



T 





vk+ 

t± 



vk- 










vk 4 



tx 

vk- 




(33) 



Due to the particle-hole symmetry, the bulk gap closing 
must occur at e = 0. Therefore, the equation should 
satisfy the following 

U 2 t\ + M 2 + (U 2 - v 2 k 2 ) 2 - M 2 [t 2 ± + 2(U 2 + v 2 k 2 )} 
= 0. (34) 

In order to have a real solution, we need 

(M 2 + U 2 ) 2 - (M 2 - U 2 )(M 2 - U 2 - t\) > 0, (35) 
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FIG. 15. (Color online) Phase diagram of bilayer graphene in the (M, U) plane at different fixed Rashba spin-orbit couplings 
Xr. Solid lines represent the phase boundary U = ±M, which separate the quantum valley-Hall region and the quantum 
anomalous-Hall region. The dashed lines are given by Eqs. (|28|) - (|29| ). separating the two different quantum anomalous Hall 
phases [see panel (a)] or the two different quantum valley-Hall phases [see panel (c)]. (a) Xr < t±. the quantum valley-Hall 
phase is characterized by C v — 4sgn(!7) (in white), while the quantum anomalous-Hall region comprises two different phases 
represented by C = 4sgn(M) (in blue/dark) and C = 2sgn(M) (in gray), (b) Xr = t±. There is no sub-phase in each region. 
The quantum valley-Hall phase and quantum anomalous-Hall phase are respectively characterized by C v = 4sgn((7) (in white) 
and C = 4sgn(M) (in blue/dark), (c) Xr > t±. The quantum anomalous-Hall region has only one phase with Chern number 
being C = 4sgn(M) (in blue/dark), while the quantum valley-Hall region includes two different phases with C v = 4sgn(£/) (in 
white) and C v = 2sgn([7) (in gray). 
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FIG. 16. (Color online) Phase diagram of bilayer graphene in the (U , Xr) plane at different fixed exchange fields M. Solid 
lines represent the phase boundary U = ±M, which separate the quantum valley-Hall region and the quantum anomalous- 
Hall region. The dashed lines represent phase boundary separating two different sub-phases in the same region, given by 
U 2 + ti - M 2 -X 2 R = 0. (a) M=0. This is just a profile along t±/t = 0.1428 in Fig. 0to compare with panels (b) and (c). 
Red/dark regions are Z2 topological insulator phase with Z% = 1 and C v = 2sgn(C7). The white region is the C v = 4sgn([/) 
quantum valley-Hall phase, (b) M < tx- Inside the region of \U\ < \M\, quantum anomalous-Hall phases emerge. For larger 
Ajj, the Chern number is C = 4, while for small Xr approaching zero, it becomes C — 2. Outside the region of \U\ < \M\, the 
valley-Chern numbers remain the same as those without exchange field as shown in panel (a), (c) M = t±. The phase boundary 
becomes linear U 2 — X 2 R = 0, which serves as a critical point changing the topology of the phase boundary, (d) M > t±. The 
essential physics remains similar to that in panel (b). 
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which gives rise to the phase transition condition 



In Fig. 1141 we provide a vivid three-dimensional (3D) 
plot of the phase diagram in the (M, Xr, U) space. For 
clarity, we do not label each phase, but will distinguish 
them in the subsequent 2D phase diagrams in detail. One 
can observe that the whole 3D space is divided by a set 
of mutual vertical planes and a uniparted hyperboloid 
determined by Eqs. (|2"8|) and (|2T))) . respectively. It is 
noteworthy that the plane of A_r = is distinct from 
other phase boundaries, i.e., the region labeled as red is 
a metallic phase. Below, we will explain the phase dia- 
gram by considering some representative regions. 

Figure [T5] exhibits the phase diagrams in the (U, M) 
plane at four different fixed Rashba spin-orbit couplings: 
(a) Xr = 0, (b) X R < t±, (c) \ R = t± and (d) Xr > t±. 
In panel (a), one observes that for small U the phase 
boundary is nearly linear to divide the metallic phase 
and quantum valley-Hall phase with valley Chern num- 
ber C v = 4sgn(i7), while for larger U the phase phase 
boundary becomes a constant. As can be seen from 
other three graphs in (b)-(d), the fundamental division 
of the parameter space into QVHI phase and quantum 
anomalous-Hall phase are separated by the solid lines 
given by U — ±M. In our calculation, the total Chern 
number is defined by C = Ck+Ck 1 - In Fig. HST b). the val- 
ley Chern number in the QVHI phase is C v = 4sgn([/) in 
the white regime, while the quantum anomalous-Hall re- 
gion comprises two different phases of matter with Chern 
numbers being C — 2sgn(M) and C — 4sgn(Af) denoted 
in gray and blue, respectively. When Xr = t±, the phase 
boundary is only determined by U = ±M , which is the 
same as we discussed in Ref. [13]. For a larger A.r as 
plotted in Fig. flST d) , the Chern number in the quantum 
anomalous-Hall phase is C = 4sgn(Af), while the QVHI 
region includes two different phases characterized by val- 
ley Chern numbers C v — 2sgn(Z7) and C v — 4sgn(£7), 
represented in gray and white. It is interesting to point 
out that at fixed M = in the gray regime, it is both a 
Z 2 TI and a C v = 2sgn([/) QVHI. 

In the above phase diagram, it is not obvious how the 
Z2 TI phase is affected by the presence of exchange field. 
In Fig. [151 we plot the phase diagram in the (U, Xr) plane 
at four fixed exchange fields: (a) M — 0; (b) M < t±; (c) 
M = t± ; and (d) M > t± . Figure USTa) shows the phase 
diagram in the absence of exchange field, which is the 
profile of t± = 0.1428i in Fig. [5] We use gray and white 
colors to denote the Z2 TI phase and conventional QVHI 
phase, respectively. When the exchange field is turned 
on, in Figs. I16f b)-(d). one finds that two different quan- 
tum anomalous Hall phases with Chern numbers C = 2, 4 
are induced when U < M at nonzero Rashba effect. It 
is noteworthy that at finite M the Z2 TI phase vanishes 
due to the time-reversal symmetry breaking, but the val- 
ley Chern number remains quantized C v — 2 sgn([/) when 
\U\ > M. From these three graphs, one observes that the 
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FIG. 17. (Color online) Bulk band gaps around valleys K and 
K' as a function of the exchange field M at fixed U = O.OSt 
and Ah = 0.15t. Here, the interlayer hopping is set to be 
t± — 0.1428t. Solid and dashed lines represent the bulk band 
gap around K and K', respectively. It can be clearly seen 
that the two bulk gaps are unequal in general, except for the 
two critical points where the bulk gaps are completely closed. 
Consistent with the phase-diagram in Fig. llSf c). it is a C v =4 
quantum valley-Hall insulator for small M before the first 
bulk gap closing; when M is located in the interval between 
two bulk gap closing points, the system enters into a C v = 2 
quantum valley-Hall phase; for even larger M exceeding the 
second critical point, it goes into a C — 4 quantum anomalous- 
Hall phase. 

phase boundary labeled by the dashed lines are governed 
by Eq. (|29|) . which reduces to the phase boundary equa- 
tion Eq. Q in the limit of M = 0. This phase boundary 
indicates a continuity with and without exchange field. 
One also observes that with increasing exchange field, the 
topology of the phase boundary in dashed line changes at 
M = t±. For zero Rashba SOC, when the exchange field 
is small, it closes the bulk band gap induced by small po- 
tential difference [see the vertical red line in panel (b)], 
driving the QVHI phase into a metallic phase; when the 
exchange field is large enough, the bulk gap from any 
potential difference is closed, giving rise to a complete 
metallic phase [see the red lines in panels (c) and (d)]. 

It is important to state that so far there are only a 
few paper a 32 ' 33 that report tunable Chern numbers in a 
quantum anomalous Hall system. From the above analy- 
sis, it is clear that our system provides another platform 
that hosts quantum anomalous Hall phases with different 
Chern numbers. The above phase diagrams are summa- 
rized concisely in Table HI which gives a complete classi- 
fication of all possible topological phases in the gated 
bilayer graphene with Rashba spin-orbit coupling and 
exchange field, shown together with the necessary con- 
ditions for a particular phase to occur. 

Another interesting feature in our system is that when 
all the three parameters U, M and Xr are nonzero, the 
bulk gaps at valleys K and K' have different responses. 
As an example, in Fig. [T7] we present the results for the 
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TABLE I. Summary of different topological phases in bilayer graphene in the presence of interlayer potential difference U, 
Rashba spin-orbit coupling An, and exchange field M. They can be divided into two categories: quantum valley-Hall insulator 
(QVHI) and quantum anomalous-Hall insulator (QAHI). Note that the Z2 Tl also belongs to the QVH1 phase. 







< |Ah| < t± 


= tx 


|Ah| > t ± 


QVHI 


M — 


C v = 4sgn(l7) 




U 2 >X\-t\: C v = 4sgn(f7) 


U 2 < A| - t 2 x : C v = 2sgn(l7) and Z 2 = 1 


< M 2 < U 2 


C v = 4sgn([/) 


M 2 < U 2 + t\ - \ 2 R : C v = 4sgn(C7) 


M 2 > U 2 + t 2 ± - A|: C v = 2 sgn((7) 


QAHI 


U 2 < M 2 


M 2 >U 2 + t\ - X%: C = 4sgn(M) 


C - 4sgn(M) 


M 2 < U 2 + t\ - \ 2 R : C = 2sgn(M) 



bulk band gaps as a function of the exchange field M 
at fixed interlayer potential difference U — 0.03t and 
Rashba spin-orbit coupling = 0.15<. It is clearly 
seen that as long as the exchange field term is turned 
on, the bulk gap amplitudes between K (solid line) and 
K' (dashed line) become unequal. However, even though 
the bulk gaps around the two valleys evolve quite dif- 
ferently, they close simultaneously. Again, the critical 
values of M at the closing points agree very well with 
the analytic expression we derived in Eqs. (j2"8")) - (f2"9"|) . 



VII. SUMMARY 

We have derived low-energy Hamiltonian descriptions 
for the Tl phase and the QVHI phase in AS-stacked 
bilayer graphene with interlayer potential U and Rashba 
spin-orbit coupling Xr. We have explored the cases when 
the bilayer graphene has an AA-stacking or has different 
Rashba spin-orbit coupling strengths in the top and bot- 
tom layers. We showed that a Z 2 Tl state can only be 
realized in the Ai3-stacked but not the AA-stacked bi- 
layer graphene. To induce a strong enough Rashba spin- 
orbit coupling in bilayer graphene, e.g., by heavy metal 
dopants or a substrate, different Rashba spin-orbit cou- 
pling strengths in the top and bottom layers Ai ^ A2 
could arise. We find that the Tl phase can be realized 
as long as AiA 2 > t\ for small interlayer potential dif- 
ference. When the time-reversal symmetry is broken by 
an exchange field M, additional topological phases can 
be induced. We find that the QVHI phase and quantum 
anomalous-Hall phase are divided by U = ±M. When 
A 7? 7^ f_i_, there exists another topological phase bound- 
ary determined by U 2 + t\ — M 2 — A|j = 0. For fixed 
\r, when \r < t±, the quantum anomalous-Hall phase 
contains two different regions characterized by the Chern 
numbers of C = 2,4sgn(M); when A,r > t±, the QVHI 
phase separates into two regions characterized by the val- 
ley Chern numbers of C v — 2,4sgn([7). Moreover, we 
find that when any two of the three parameters (i.e. in- 
terlayer potential difference, Rashba spin-orbit coupling, 
and the exchange field) are considered, the bulk band 
gaps at K and K' are equal. However, if all three terms 



are present, the bulk gaps at K and K' become different 
except at the topological phase transition points. It is 
noteworthy that in multilayer graphene, a bulk band gap 
opens in the presence of an external electric field. This 
makes multilayer graphene a good candidate to explore 
more interesting topological phasesi« 
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A. Appendix 

In the following discussion we set Xr and t± as fixed, 
and allow U to slightly deviate from the phase transition 
point Uo, i.e., U = Uo ± A, where Uq = \J A^ — and ± 
correspond respectively to the quantum valley-Hall and 
topological insulator phases. 

As stated in the main text, when the topologi- 
cal phase transition occurs, the bulk band gap closes 
at the exact Dirac points K/K'. At the crit- 
ical Uo, the system Hamiltonian on the basis of 
{B u , A 2 f, Ai t , A 2i,B 2 f, B 2 i, B n , A n } can be expressed 
as 



H(U ) = 





























-Uo 
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-Uo 
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t± 




Uo 



Through a direct diagonalization, the eigenenergies are 
obtained as e = 0, 0, ±Uo, ±£1, ±e 2 , where S\ — —Uo/2 + 
y/8X 2 R + U'^/2 and e 2 = -U /2 - y/8X 2 R + U^/2. The 
former four correspond to the low-energy part, while the 
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latter four correspond to the high-energy part. Based on 
our analysis of the Berry curvature from the tight-binding 
model, the valley-Chern number arises only from the low- 
energy bands and the high-energy bands contribute zero. 
For the diagonal block Hamiltonian of h\(Uo) = 

Uo t±_ 

its unitary transformation matrix is 



which leads to 







-Uo 



Vu V12 V13 

V21 V 2 2 "23 

V31 «32 V33 

t± 



l/a 2 



l/a 3 



V2\r 

i »Ar(C/ -£i) »A R (L r -e 2 ) 

y/2 a 2 t ± (U +ei) a 3 t ± (U +e 2 ) 
-Up Ei— Up £2-Up 
\f2\n ct 2 t± a 3 t± 

which can diagonalize h\(Uo) to be 



h' 1 (U )=V\h 1 (U )V 



ct2 and CK3 are respectively 




£1 
e 2 



OL2 



0L3 




Ar(?7o — £ i)i9 , r £ i — ^Oi 
] 2 + [- 



(Uo + £i)*_l 



t± 



/ 1 + jAr(J7o_£2) j 2 + ^2 - Up ^ 



(Uo + e 2 )t± 



t± 



(37) 



(38) 



For the other diagonal block Hamiltonian of h 2 (Uo) 



-Uo 





t± 







Uo 


i^R 


, its unitary transformation matrix 


t± 


-i^R 


Uo 





is written as 



V, = 



V\l V12 V\3 
-V21 -V 2 2 -V23 
-V31 -V32 -V33 



(39) 



h' 2 (Uo)=Vlh 2 (Uo)V2 





-e x 
-E 2 



In order to arrange the eigencnergies to be low and 
high energy parts, the basis should be reordered to 
be: {A n ,B 2 i,B n , A 2 f, A 2 ^, B 2 f, Bi T , A n }. The corre- 
sponding full unitary transformation matrix becomes: 






vu 
V21 
V31 

















-V21 
-V31 
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V12 V13 

v 2 2 v 2 3 

U32 W33 















V12 
-V21 









V\3 

-V23 



-V32 -V33 



When the interlayer potential difference is slightly de- 
viated from Uo, the Hamiltonian is written as: 



H(U) = H(U ) + H A + H k , 



(40) 



where 



H A = Diag{A, -A, A, -A, -A, -A, A, A}. (41) 



and 



Hk 
























vk. 














vk- 





























vk- 




















vk- 











vk + 





























vk + 




















vk + 

















vk- 
























By performing a unitary transformation, the Hamilto- 
nian of H(U) becomes 
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H'(U)= V t i7(/7)V 

(2ufj - 1)A 7i«fc_ v 3 ivk + 



-livk + 


(l-2«?i)A 


-w 3 iwfc_ 








-u 3 iwfc + 


A + C/o 





vzivk- 








-A-t/o 


2u u v 12 A 


7 2 wfc_ 







2v lx vi 3 A 


7 3 vfc_ 





V 33 vfc + 


~l2vk + 


-2wiit>i 2 A 


-v 32 vk_ 





-J3vk + 


-2uii«i 3 A 


-v 33 vk- 






Hp T 



2W11V12A 2viit>i 3 A 72«fc- "f 3 vk- 

-7 2 vfc+ -7 3 wfc + -2vn«i 2 A -2wnVi 3 A 

-w 32 i;fc + -w 33 i;fc + 

v 3 2vk- v 33 vk~ 

74 * A + Ei 75A — 2w 12 «22Vfc_ 78^^- 

75A 7 6 A + e 2 78«fc- -2wi 3 w 23 i;A:_ 

2wi 2 w 22 wA; + -^ 8 vk + -74 * A - e\ 77A 

-j 8 vk + 2ui3f; 2 3u/c + 7 7 A -76A - e 2 



(42) 



where 71 = —2vuv 2 i = -it±/\ R , j 2 = -(W11U22 + 

Wl 2 f2l), 73 = -(Ull«23 + Wl3«2l), 74 = U?2 + U 22 ~ u 32 ) 
75 = V 12 V 13 + V22V23 ~ W32W33, 76 = «13 + «23 ~ w 33' 

77 = W32W33-W22W23-W12W13, and 7 8 = -(vi 2 t>23+vi3^22)- 
Since both A and fc + / h- are extremely small, the effec- 
tive Hamiltonian can be simplified to be 



Explicitly, Hp can be written as: 



Hp = 



itj_vk+ 



— Uovk- 
V2\r 



-it±vk- 






A + C/q 







-A-Uo 



(44) 



For small A and k, the higher energy block can be 
expressed as: 



H 



eft 



Hp - TH7, l TK 



(43) 



Hn = 



£1 

£ 2 










-El 









-£2 



(45) 
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